Mathematical Biosciences

Immune cell engagers (ICE) such as bispecific antibodies (bsAbs), within an immunological synapse, bind and link CD3 on a T cell to a target antigen (TAA) on a cancer cell, forming a trimer (CD3:bsAb:TAA complex). With sufficient trimer numbers within the synapse, the T cell can become activated and promote cancer cell killing. Elranatamab, a CD3-bispecific antibody for multiple myeloma, has received FDA and EMA filing acceptance (August 2023 and December 2023, respectively) adding to a growing list of bsAbs that are treating patients. In the drug development stages of ICE bsAbs, mechanistic modeling approaches are often used to attain a greater quantitative understanding of the modality, preclinically, and provide human pharmacokinetic and efficacious dose predictions to aide in Phase 1 trial design. To date, the majority of ordinary differential equation (ODE) trimer models treat the tumor compartment as well-mixed and trimer formation is governed by a bulk population reaction not accounting for individual synapses. This lack of discrimination can lead to imprecise analysis when analyzing results across E:T ratios using metrics like trimers per T cell or trimers per target cell. To this end we developed an ODE trimer model based on single-synapse complexes (one target cell/one immune cell) with 2D cross-linking trimer formation. We show computationally that the number of trimers per synapse is invariant to the value of the E:T ratio for a given free bsAb concentration, a property that cannot be captured by non-synapse models. A simple demonstration of this discrepancy using the well-known Betts trimer model is presented. We then apply the Betts trimer model coupled to a tumor growth inhibition (TGI) module to show that our synapse-based trimer model is easy to substitute in to model TGI, including the addition of a trimer-per-synapse activation threshold function for cell killing. Overall, our model attempts to balance mechanistic fidelity while limiting the complexity of the model.

Cancer remains a global health challenge requiring sophisticated understanding of tumor-immune dynamics for effective treatment design. Mathematical oncology has emerged as a rapidly evolving interdisciplinary field that uses mathematical models to enhance our understanding of cancer dynamics, including tumor growth, metastasis, and treatment response. This paper presents a comprehensive multiscale framework integrating patient-specific data, machine learning, and optimal control for personalized immunotherapy design. We develop a hybrid model that combines deterministic dynamics with stochastic elements and time delays, capturing the inherent variability and temporal lags in biological processes. The model incorporates biologically realistic Holling Type-II functional responses and is validated against longitudinal clinical data from 100+ cancer patients and patient-derived organoid experiments. Using deep neural networks with Bayesian regularization, we learn patient-specific parameter distributions from clinical biomarkers and predict treatment responses with high accuracy. Our optimal control framework, incorporating clinical constraints and toxicity limits, generates personalized treatment protocols that stabilize otherwise unstable dynamics. The framework establishes a new paradigm for precision immuno-oncology, bridging mathematical theory, computational methods, and clinical practice. Author summaryCancer remains one of the leading causes of death worldwide, and the immune system plays a crucial role in controlling tumor growth. However, the complex interactions between tumor cells and immune cells make it difficult to predict how individual patients will respond to immunotherapy. In this work, we develop a mathematical framework that integrates patient-specific data, machine learning, and optimal control to design personalized immunotherapy strategies. Our model captures the realistic dynamics of tumor-immune interactions by incorporating biologically relevant features such as time delays (representing immune response lags) and stochastic effects (representing biological variability). Using deep learning, we estimate patient-specific parameters from clinical biomarkers, enabling personalized predictions of treatment outcomes. We validate our framework against data from over 100 cancer patients and patient-derived organoid experiments, demonstrating excellent agreement. Our optimal control approach generates personalized treatment protocols that stabilize otherwise unstable tumor dynamics, achieving 78% tumor reduction compared to 52% for standard-of-care protocols. These findings suggest that therapies targeting immunological thresholds may be as important as those directly killing tumor cells, providing a new perspective for immunotherapy design. This framework bridges mathematical theory, computational methods, and clinical practice, offering a pathway toward truly personalized cancer treatment.

Auxins are a family of plant hormones involved in various processes across plant tissues and species. The Nuclear Auxin Pathway (NAP) consists of interacting transcription factors (ARFs) and repressors (Aux/IAAs), which govern an individual cells response to changes in auxin concentration. These components are present in all land plants, and many species possess multiple copies of each signalling component. We present a general framework for ODE-based models of NAP submodules with the flexibility to model the promotion and repression of target genes by any combination of transcriptional regulators. We analyse published data and show that auxin treatment in Arabidopsis thaliana roots triggers a range of characteristically distinct temporal response profiles--for both target genes and the signalling components themselves. Using our modelling framework, we recapitulate aspects of this behaviour by presenting examples of real and theoretical NAP subnetworks, and by analysing the effect that these network dynamics have on auxin-mediated transcriptional responses. This work demonstrates the utility of our modelling framework as a general-purpose tool for understanding the function of certain protein-protein and protein-DNA interactions through their effects on the NAP. This exploration of the rich dynamics of more complex signalling pathways promises to advance our understanding of the NAP.

1.Partial Rank Correlation Coefficient (PRCC), usually performed following Latin Hyper-cube Sampling (LHS), is a global sensitivity analysis that quantifies the monotonic relationship between model parameters and the desired output. To carry out this analysis, a range of acceptable parameter values must be known or estimated. However, within a biological context, approximating these values may be difficult. Parameter values and ranges can be taken from different organisms or systems or be estimated to produce qualitative phenomena in the model. Using a mathematical model of the epithelial mesenchymal transition (EMT) as a test case, this work examines how the parameter ranges chosen prior to analysis can influence LHS-PRCC results and shape subsequent analysis interpretations. Previous LHS-PRCC analysis of this model restricted parameters to {+/-}10% of their original value, which limits the scope and interpretability of parameter influence. Such a small range assumes, in the biological sense, that parameters are well-measured with little variability. Here, this work extends the previous analysis and explores several parameter ranges ({+/-}25%, {+/-}50% of the original value). This work also tests whether, within the {+/-}10%, {+/-}25% and {+/-}50% parameter ranges, the bistable switch present in the original model are maintained. Ultimately, this work showcases how a choice made prior to analysis, such as the accepted parameter ranges for biological rates and values in complex dynamical systems can influence sensitivity analysis results and interpretability. Additionally, these choices can have hidden consequences, such as the loss of phenomenological behavior. Thus, explicit prior knowledge about the appropriate parameter values is needed before using analysis to guide future experiments and model development.

In a population model for an infectious disease, we consider the early stochastic dynamics of an emergent 'mutant' strain, appearing and spreading during an epidemic of another 'wildtype' strain. The mutant may not reach establishment in the host population. The time at which the mutant first appears determines its probability of establishment. We calculate this establishment probability with two methods. The first method assumes a classical branching process, with a constant transmission rate. The second method reflects the changing size of the pool of susceptible hosts, due to the dynamics of the wildtype. We find that susceptible depletion can substantially impact the establishment probability. We explore the consequences of this stochastic establishment on the "escape pressure" acting on a pathogen to produce immune escape variants. We find that the overall escape pressure rate depends strongly on the appearance time of the mutant, especially if the establishment probability is itself shaped by the continued spread of the wildtype. In most scenarios, the escape pressure rate (and thus, the risk of new escape variants) peaks slightly earlier than the prevalence of the wildtype strain. Integrating the escape pressure over time, we obtain the cumulative escape pressure generated by the wildtype epidemic. The relationship between the escape pressure and the vaccination coverage depends on the cross-immunity, due to susceptible depletion. For example, with intermediate cross-immunity, the risk of immune escape may be lowest at intermediate vaccination coverages. Thus, these results raise important considerations for vaccination strategies in response to novel outbreaks.

Diffuse gliomas are deadly because the individual tumor cells invade - they travel far from the imageable mass, so it is impossible to remove the tumor completely. On the cellular level, glioma cells seem to be in either a "go" state (in which they do not divide) or a "grow" state (in which they do not migrate). We investigate what this tiny choice has to say about the large-scale speed of the invasion front and whether the implication is sufficiently strong to rule out the classical description of the Fisher-Kolmogorov-Petrovsky-Piskunov (Fisher-KPP) type, in which a single phenotype migrates and proliferates. We derive a two-phenotype reaction-diffusion model with density-dependent switching, and we prove the cooperative (quasi-monotone) structure and the associated comparison principle and study travelling-wave solutions of the model. A leading-edge linearization gives minimal front speed as minimizer of an explicit dispersion relation, and direct simulation verifies the predicted speed. In the experimentally relevant fast switching limit, we find a closed-form expression for the speed, that is, we obtain an effective Fisher-KPP equation with rescaled diffusivity and growth rate, with the fractions of the phenotypes. The "go-or-grow" (GoG) front can move at a maximum speed of half the Fisher speed for the same single-cell motility $D$ and proliferation rate $r$, which occurs only when the cells divide their time equally between the two phenotypes. This bound is directly testable: measurement of the front speed, plus independent determination of $D$ and $r$, discriminates the two hypotheses, and in the GoG case, yields recovery of the phenotype balance. We then extend the result to anisotropic (DTI-informed) invasion along white-matter tracts and discuss implications for understanding clinical measurements of growth rate.

Dengue transmission is shaped by multiple viral serotypes, temporary cross-immunity (TCI), antibody-dependent enhancement (ADE), and repeated exposure in endemic populations. Classical multi-strain models usually assume lifelong protection against reinfection with the same serotype. However, recent evidence suggests that homologous dengue reinfections, although rare, can occur. Their population-level consequences remain poorly understood. We extend a two-infection, two-strain dengue model with TCI and ADE-mediated transmission differences to include homologous reinfections. Homologous reinfection is represented by two exploratory parameters: relative susceptibility to reinfection with the same serotype and relative infectiousness during homologous reinfection. Using equilibrium analysis, bifurcation diagrams, simulations, and phase-space projections, we examine how these parameters affect dengue dynamics and interact with TCI duration and seasonal forcing under intermediate and long TCI durations, with and without seasonality. The extended model shows that qualitative dynamics characteristic of endemic dengue transmission are reproduced mainly when susceptibility to homologous reinfection is low, so that homologous reinfections remain rare but dynamically influential. Longer TCI broadens regions of complex oscillatory dynamics, while seasonality shifts the bifurcation structure and makes torus bifurcations a central route to complex behavior. Although backward bifurcation can occur when homologous susceptibility exceeds the biologically meaningful range, this result should be interpreted as a mathematical mechanism rather than a realistic dengue scenario. These results indicate that rare homologous reinfection pathways can influence long-term dengue dynamics when interacting with immune history, TCI, ADE-mediated transmission differences, and seasonal variation. Incorporating such pathways may improve understanding of recurrent outbreaks and irregular incidence patterns in highly exposed populations.

Asymptotic standard errors for the parameters of a nonlinear mixed-effects model fitted by first-order conditional estimation (FOCE) or FOCE with interaction (FOCEI) require the observed (Fisher) information -- the negative second derivative of the population objective at the optimum. The gradient of this objective can be computed exactly from sensitivity equations, but the observed information is conventionally still formed by finite differencing, which is less accurate and step-size dependent. Our objectives are to (i) derive the FOCE and FOCEI observed information in closed form within the same sensitivity-equation framework, and (ii) quantify the precision this recovers. Writing the objective as a data term plus the log-determinant of the first-order inner Hessian, the population Hessian splits so that the data term reuses the second-order sensitivities already needed for the gradient, whereas the log-determinant term requires third-order sensitivity equations -- confining the third-order dependence to a single term, where it enters in exactly two places. A finite-difference error analysis shows the differenced Hessian attains an accuracy no better than the square root of the objectives evaluation accuracy, whereas the analytic form is limited only by the sensitivity and differential-equation solutions, with no step size to tune. We illustrate this on a one-compartment oral model with first-order absorption fitted to warfarin data, where the differenced standard error is usable only over a narrow band of step sizes while the analytic value carries none. Implemented in the open-source R package nlmixr2, the method makes exact, reproducible standard errors routine, supporting more dependable confidence intervals, identifiability assessment, and uncertainty propagation.

Malaria and diabetes represent two globally significant metabolic disorders whose co-occurrence leads to complex, poorly understood pathophysiological interactions. Plasmodium infection disrupts glucose homeostasis through parasite-driven glucose consumption, inflammatory cytokine production, and pancreatic /{beta}-cell dysfunction, while diabetes impairs host immunity and increases malaria susceptibility. To date, no mathematical framework has captured the bidirectional coupling between these systems. Here we extend the insulin-glucose-glucagon (IGG) model of Dalton et al.\ (2026) by introducing a fourth state variable representing parasite load, incorporating malaria-induced insulin suppression, parasite-driven glucose consumption, inflammatory gluconeogenesis, bidirectional glucagon dysregulation, and insulin-dependent immune enhancement of parasite clearance. We establish positivity, boundedness, existence and uniqueness of steady states, local stability via Routh-Hurwitz criteria, global stability via Lyapunov functions, and sensitivity analysis of parameters driving hypoglycemia risk. Numerical simulations characterise the model across healthy, diabetic, and co-infected states. They show that parasite-driven glucose consumption and inflammatory gluconeogenesis act antagonistically on circulating glucose, that insulin-enhanced immunity lowers peak parasitemia through a saturating clearance term, and that increasing the half-life of exogenous insulin raises hypoglycemia risk in all host states. These mechanisms provide testable hypotheses for the clinical management of malaria-diabetes patients and identify potential therapeutic targets (TNF- blockade, glucagon analogues) for mitigating co-infection morbidity.

The global potential of a chemical reaction network has many applications and is closely related to the stochastic detailed balance. However, many fundamental questions concerning stochastic detailed balance remain unresolved, such as whether it depends on the system volume and how to construct new systems that satisfy it. In this paper, we show that stochastic detailed balance may depend on the system volume. We therefore introduce four types of stochastic detailed balance according to their dependence on volume and rate constants, and systematically investigate the relationships among them. Our results distinguish detailed balance arising from particular choices of volume and parameters from that enforced by network structure, and identify conditions under which detailed balance at one volume extends to all volumes. We further obtain a class of networks satisfying stochastic detailed balance for every volume and every positive choice of rate constants, and construct new systems whose global potentials exhibit double-well structures.

Background: Tuberculosis remains a persistent infectious disease whose control is complicated by latent infection, delayed treatment, incomplete recovery, reinfection, and continuing transmission from infectious individuals. Although treatment is central to tuberculosis management, it is frequently represented only as a transition parameter in mathematical models rather than as a separate epidemiological state. In this study, treatment was therefore incorporated explicitly as an independent compartment so that its influence on transmission, recovery, disease-induced mortality, and long-term disease persistence could be evaluated. Methods: A deterministic nonlinear compartmental model was formulated by dividing the total population into susceptible, exposed, actively infected, treated, and recovered classes. Reinfection of recovered individuals, progression from latent infection to active disease, movement of infectious individuals into treatment, treatment-associated recovery, natural mortality, and disease-induced mortality were included. Positivity and boundedness of the solutions were examined to establish biological validity. The basic reproduction number, R0, was derived through the next-generation matrix approach. Disease-free and endemic equilibria were determined, and their local and conditional global stability properties were investigated using Jacobian analysis, the Routh-Hurwitz criterion, center manifold theory, Lyapunov functions, and LaSalles invariance principle. Normalized sensitivity indices, Latin hypercube sampling, partial rank correlation coefficients, and numerical simulations were also applied. Results: The disease-free equilibrium was shown to be locally asymptotically stable when ,R0<1 whereas sustained transmission and a unique endemic equilibrium were associated with R0>1. Under the stated reduced-model assumptions, stability of the endemic equilibrium was established. Transmission-related parameters were identified as the strongest positive contributors to disease persistence. In contrast, treatment and recovery parameters were found to reduce the reproduction number and infectious burden. Numerical simulations indicated that stronger treatment implementation and reduced transmission opportunities produced substantial reductions in active tuberculosis cases. Conclusion: Treatment was shown to function as both a clinical pathway and an epidemiological control mechanism. The proposed framework may support the design of treatment-centered strategies for reducing tuberculosis prevalence and preventing long-term endemic persistence.

We develop a mathematical model of Rift Valley Fever integrating mosquito vectors, ruminants, and humans, based on an SEIR-type structure with vertical transmission in vectors. Local data from the Sudanian and especially the Sahelian zones are used to capture the impact of climatic variations on mosquito population dynamics. The mathematical analysis establishes the models positivity, determines the basic reproduction number R0, and demonstrates the local and global stability of the disease-free equilibrium. Sensitivity analysis (PRCC) highlights the most influential parameters, while the stochastic approach using a continuous-time Markov chain confirms the major role of seasonal rainfall. Numerical simulations reveal a peak in animal and human infections around the 9th month, correlating with periods of heavy rainfall. This model provides a relevant tool for surveillance and prevention within a "One Health" approach in Chad.

Immune checkpoint blockade (ICB) therapy has revolutionized cancer treatment, though it is still effective in only about 20-40% of patients. To increase the efficacy of ICB therapy, combinations of ICB antibodies such as anti-PD1 and anticoagulants (e.g., thrombin inhibitors) have been studied in preclinical mouse models and clinical trials. We developed a Bliss-type analysis to quantify the synergy between anti-PD1 and dabigatran etexilate using published tumor growth data from a mouse model. We then developed a minimal mechanistic computational model to quantitatively study the synergy between anti-PD1 and the thrombin inhibitor dabigatran etexilate in a published study (Metelli et al.) of a preclinical mouse model of colon cancer. Our model included tumor cells, CD8+ T cells, and the pleiotropic cytokine TGF{beta}, whose production in platelets is influenced by thrombin, and described a potential mechanism of interplay among these components in the tumor microenvironment. We performed nonlinear mixed-effects modeling to capture mouse-to-mouse variation in tumor growth under different treatment conditions. The estimated parameter values pointed to several underlying mechanisms of synergy, including increased expansion of CD8+ T cells in the tumor microenvironment in the presence of anti-PD1 and dabigatran etexilate. These predictions can be further validated in future experiments. Thus, the combination of longitudinal tumor growth data, mechanistic population dynamics modeling, and nonlinear mixed-effects modeling can be used to study synergy between ICB and other drugs in controlling tumor growth.

Oscillations and oscillatory synchronization are pervasive in motor circuits, where their role in rhythm generation and entrainment is well established but their role in feedback control of movement remains unclear. Here I show analytically that two oscillators of any type, coupled through a delayed interaction that is an odd function of their phase difference, necessarily implement a proportional-derivative (PD) control law in the near-synchrony limit. The proportional gain follows from the slope of the coupling function and the derivative gain is set by the coupling delay, so that PD control emerges with no additional machinery. Simulations confirm that such oscillators reproduce ideal PD step responses near synchrony and that control quality degrades systematically away from it. This establishes a direct, model-independent bridge between oscillatory synchronization and feedback control, and suggests concrete experimental signatures for candidate systems.

In the face of an epidemic where a population behaviour both influences disease transmission and reacts to it, social processes can generate norms to support socially beneficial behaviour. Most mathematical models of coupled behaviour disease dynamics treat norms as pre-existing rather than explaining how they are maintained. Here, we investigate whether altruistic punishment can sustain a social distancing norm when individuals may defect, cooperate without punishing, or cooperate while paying a cost to punish defectors. We couple a transmission model to an imitation model for these three strategies. Disease prevalence affects behavioural payoffs, while the behavioural composition modifies transmission. We also compare this coupled system with a control where behavioural decisions respond to a fixed prevalence. We find a wide parameter regime corresponding to the establishment of an injunctive social norm in support of social distancing, where the punisher strategy is widespread. Persistence may occur through stable states where punishers or dominant. Disease behaviour feedback can also create oscillations (where the three strategies succeed one another in response to epidemic waves) or tipping points (sharp transitions between all-defector and cooperative states). These effects do not occur in the uncoupled model, although there are still broad parameter regimes where a social norm persists. Our findings show that costly peer punishment can support persistence of social norms that mitigate disease transmission. More broadly, endogenous epidemic feedback can qualitatively change the conditions under which cooperation and punishment are sustained, producing tipping points and long-term behavioural epidemiological cycles that fixed-payoff models cannot capture.

Avian Influenza viruses (AIVs) infect a broad host range despite having a natural reservoir in wild aquatic birds. Whilst most strains stay within their host species, some break the species barrier through genetic adaptations. We are most concerned about zoonotic cases, where a human becomes infected. Despite these events being rare, they are associated with high mortality and introduce the risk of onward human-to-human transmission of AIV. As a novel pathogen within the human population, this could have pandemic potential. Using genetic composition features for 8 AIV proteins drawn from viral sequence data, we employ machine-learning algorithms to classify AIV cases as zoonotic or not. These genetic features encode host 'signatures' which can indicate zoonosis and include frequency measures such as dipeptide composition and amino acid physiochemical properties. We consistently find XGBoost to outperform all other algorithms. We optimise parameters for ten classification models: one for each of the 8 proteins and two combined models. Following this, we show that a multi-model approach gives the best performing prediction for AIV zoonosis. We have identified all 8 proteins as having a role in predicting zoonotic transmission. Of particular importance is the PB2 and HA proteins, with specific amino acid physiochemical properties such as charge, secondary structure and hydrophobicity amongst the most indicative features in our combined models. Our alignment-free computational study can identify AIV cases still within avian hosts which are genetically closest to zoonotic AIV cases, thereby identifying the cases most likely to cross the species barrier. In a resource limited environment, our model could be used to quickly identify high priority cases for further investigation.

The basic reproduction number R0 confounds pathogen biology with adaptive human contact behavior. Earlier epidemiological--economic theory predicted a forward-looking behavioral contact response but could not test it in the absence of appropriate behavioral data. Using directly measured mobility as an observable proxy for contact, we (i) estimate the behavioral response function directly from data; (ii) show that the biology/behavior decomposition and hence the behavioral correction to R0 is not identified from an epidemic trajectory, the apparent constant-contact R0 being one endpoint of an observational-equivalence class that fits the factual curve identically yet diverges under counterfactual; and (iii) characterize that divergence ("what R0 deletes") as state-dependent, unimodal in counterfactual severity and vanishing when behavior saturates. We then show that, across US jurisdictions, the correction is empirically bounded because risk-responsiveness and behavioral non-saturation are confounded (r=-0.57, n=51): where behavior could compensate, it was already maximal, and where it was not maximal it did not respond. What R0 deletes is thus real and structurally characterizable yet empirically modest here, for reasons the framework itself supplies.

AO_SCPLOWBSTRACTC_SCPLOWCell fate decisions are driven by gene regulatory networks (GRNs). While the mutually inhibitory toggle switch effectively models binary fate decisions, fully connected inhibitory networks with more than two nodes fail to capture multi-fate decisions due to the low prevalence of "single high states", where only a single master regulator is highly expressed. The goal of this study is to find network structures that support all single high states. We find that the only network that attains the highest possible prevalence of all single high states within the set of monotone Boolean (MB) models is completely disconnected. Since biological networks typically require connectivity, we investigate network structures that support equipotency, where all single high states have equal prevalence within MB models. Finally, we characterize the networks that support multistability between all single high states, finding that it is possible only in networks in which each node either has self-activations or is inhibited by every other network node. Our findings provide a theoretical framework for understanding the network design principles that can support simultaneous differentiation into multiple distinct cell types.

Diaphorina citri is the primary vector of Huanglongbing (HLB), a devastating disease that affects global citrus production. Effective biological control using the parasitoid Tamarixia radiata represents a sustainable alternative to chemical insecticides, but its efficacy depends heavily on environmental variables and resource availability. In this work, we develop a four-population mathematical model to analyze the dynamics of D. citri in the presence of T. radiata, a known parasitoid. Our model incorporates the cyclical and periodic behavior of citrus flushing (new shoots), providing a realistic representation of resource-limited dynamics as observed in field conditions. Through comparative simulations of three scenarios, without parasitoid, a single initial introduction, and periodic augmentative releases, we characterize and compare the impact of T. radiata as a biological control agent. Our results show that periodic releases maintain pest suppression, whereas a single introduction only delays pest recovery. By aligning theoretical modeling with ecological reality, this framework supports the role of T. radiata in pest suppression and provides a baseline for future work on optimal and cost-effective release strategies.

Time compaction theory is a general framework explaining how a brain can efficiently deal with dynamic situations occurring in, e.g., sports games. It involves a geometric representation of the time dimension, which enables effective learning and strategic action planning. The theory has recently received experimental support in humans. However, its current computational model has an important limitation: it does not account for deliberate waiting and speed modulation, behaviors ubiquitous in natural environments. This work substantially extends the original model formulation by a dimensional lifting of an n-D workspace into (n + 1)-D mental space, where time remains geometrically embedded. The proposed biologically inspired computational model can generate adaptive behavior across increasingly complex situations, from navigation in everyday social environments to competitive sports. Furthermore, by actively conditioning the expected responses of other agents and stabilizing future predictions, we introduce the concept of uncertainty points in sequences of generalized cognitive maps to support the generation of adaptive strategies in interactive environments, where future prediction has a limited time horizon. Thus, we provide a mechanism for chaining short-term solutions into long-term strategies, which is illustrated by simulating the behavior of a player in a real football game. Author summaryHumans often anticipate future interactions in dynamic environments. Many behaviors, such as avoiding other pedestrians, letting someone pass through a narrow corridor, or reproducing the kind of dribbling maneuvers performed by elite football players, require deciding not only where to move but also when to move. Existing theories suggest that the brain simplifies such situations by representing future interactions as static spatial maps, making them easier to learn and recall. However, current computational models cannot naturally account for common behaviors such as waiting, slowing down, or modulating speed. Here we show that these behaviors readily emerge if the model space is extended by an additional virtual coordinate that encodes accumulated waiting rather than physical time. The proposed model simultaneously admits a wide variety of behaviors, including speed modulation, multigoal decisions, and compound actions, while preserving the principles of time compaction. We illustrate the model in everyday situations and by reproducing two real football plays, comparing the observed behaviors with model simulations. Our results suggest computational principles through which the human brain may efficiently represent, memorize, and exploit dynamic situations.